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Abstract 

Forty years ago, Robert May questioned a central belief in ecology by proving that sufficiently large or complex 
ecological networks have probability of persisting close to zero. To prove this point, he analyzed large networks in 
which species interact at random. However, in natural systems pairs of species have well-defined interactions (e.g., 
predator-prey, mutualistic or competitive). Here we extend May's results to these relationships and find remarkable 
differences between predator-prey interactions, which increase stability, and mutualistic and competitive, which are 
destabilizing. We provide analytic stability criteria for all cases. These results have broad applicability in ecology. 
For example, we show that, surprisingly, the probability of stability for predator-prey networks is decreased when we 
impose realistic food web structure or we introduce a large preponderance of weak interactions. Similarly, stabihty is 
negatively impacted by nestedness in bipartite mutualistic networks. 



May showed mathematically that large and complex ecosystems are inherently unstable fV^. This contribution 
has been one of the main drivers of theoretical ecology ever since |[3]|4]|5l, as it clashed with the prevailing belief of 
ecologists that large, highly complex ecosystems (such as those observed empirically) were more stable than simpler 
ones (found in extreme environments and disturbed ecosystems) ||6l|7]. 

May's theorem deals with a particular type of community matrix |I8] [B 121 M, of size 5 x 5* (S" is the number 
of species in the system). The matrix M describes the effect a species j (column) has on species i (row) around the 
equilibrium point of an unspecified dynamical system describing the density of the species through time. 

The diagonal coefficients of M are all —1, while the off -diagonal coefficients are drawn from a normal distribution 
A^(0, cr^) with probability C and are otherwise. For large S, May proved that the probability of stability is close 
to whenever the "complexity" K = ay/SC > 1. Local stability measures the tendency of the system to return to 
equilibrium after small perturbations. In unstable systems, even infinitesimal perturbations will make the system move 
away from the equilibrium state, potentially resulting in the loss of species. Thus, it should be extremely improbable 
to observe rich (large S) or highly connected (large C) ecosystems persisting through time. Mathematically, an 
equilibrium point is stable if all the eigenvalues of the corresponding community matrix have negative real part. 

The networks described by these matrices have random structure: each pair of species interacts with a given 
probability. However, this randomness translates, for large S, into fixed interaction frequencies, so that when we are 
constructing the matrices above, we are following a precise mixture of interaction types. We can classify interactions 
types according to the signs of the ordered pair {Mij , Mji) (effect of j on i and vice versa), and compute their expected 
frequencies in a large random matrix. For the pair (0, 0): non-interacting, the expected frequency is (1 — C)^; (+, 0) 
or (0, +): commensalism, C(l — C); (— , 0) or (0, — ): amensalism, C(l — C); (— , — ): competition, (+, +); 

mutualism, (+, — ) or (— , +): predator-prey, 

Here we show how the criterion for stability changes when we impose a specific type of interaction between 
species. We start with predator-prey matrices, which are like random matrices but with the constraint that if AUj > 0, 
then Mji < 0: the interaction is beneficial for one species and detrimental for the other Numerical simulations 
showed that these matrices are more stable than random ones f9l. This is confirmed by our results, as we find that the 
stability criterion becomes K < 7r/(7r — 2) sa 2.75. Thus, stable predator-prey systems can be much larger and more 
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complex than random ones. For example, for a — 0.5, C = 0.1, the stability criterion is violated for S" = 40 in the 
random case, but S = 304 for the predator-prey case. 

It should be noted that the mean of the off-diagonal coefficients Mij is in both cases. In fact, both random 
and predator-prey matrices will have (on average) the same number of positive and negative coefficients, and with 
the same magnitude. The only difference between the matrices is that in the predator-prey case, coefficients are 
arranged in pairs, such that one is negative and the other is positive. However, this arrangement modifies the expected 
interaction strength product for two interacting species MijMji (i.e., the expectation is taken over all the pairs in 
which Alij or Mji are ^ 0), which is in the random case, but — (T^2/7r in the predator-prey case. Thus, the 
difference in stability arises from having negative MijMji. This is confirmed by showing that the mean eigenvalue 
A = — 1 in both cases (the trace being —S), while the variance ifTOl (for large S) is Var(A) = in the random case 
and Var(A) = {S — l)CMijMji in the predator-prey case. Note that the variance can be negative given that the 
eigenvalues can be complex conjugate. Having negative variance means that the variance of the imaginary part of the 
eigenvalues is larger than that of the real part. If the stability is driven by having negative MijMji, reversing its sign 
should decrease stability. 

The matrices yielding the opposite sign for MijMji, compared to the predator-prey case, are those in which 
pairs of species interact as mutualists or competitors, and for each pair the interaction type is assigned at random. 
In these matrices, we still have Mij — 0, but now MijMji = a^lj-K. Accordingly, the stability criterion becomes 
K < 7r/(7r + 2) w 0.61: this mixture of competition and mutualism leads to a large decrease in stability. For a = 0.5, 
C = 0.1 the criterion is violated for S = 15. 

We derived these criteria for stability in the following way. Consider a random, S x S matrix. A, whose elements 
are all Gaussian with mean Aij ~ 0, mean square value A?^ — 1/S, and mean interaction strength product AijAji = 
t/S. For S ^ CO, the eigenvalues of A, \ ~ x + iy, are uniformly distributed in the ellipse (x/a)'^ + {y/b)'^ < 1, 
with a = (1 r) and 5 = (1 - r) ifTll . 

To obtain the community matrices we are interested in, rescale the matrix A: aVSA — M. Thus, the elements 
Mij have the following properties: Mij — 0, Mfj = a^, and MijMji — to^ . For large S, the eigenvalues of M are 
approximately uniformly distributed in an ellipse with a ~ a\fS{\ + t) and h = (t\/S'(1 — r). 

The value of r can be derived for all the types of matrices illustrated above. In the random case, MijMji = 
T(T^ = 0, and thus a = b = aV~S (i.e., the eigenvalues are distributed in a circle). For the predator-prey case, 
the expectation for the product MijMji is — (T^2/7r: the expectation for the product of two independent, identically 
distributed half-normal random variables, with a negative sign accounting for the opposite signs of the coefficients. 
Thus, MijMji = Ta^ = -a^lj-K, leading to r = -2/7r, a = ctV^(1 - 2/7r) and b = a^fS(\ + 2/7r). Similarly, for 
the mixture of competition and mutualism we have MijMji = to^ ~ a^lj-K and thus t — 2/11. 

These results hold for completely connected matrices, with ellipses centered at (0, 0). Setting the diagonal coeffi- 
cients to — rf centers the ellipses at (— rf, 0) (for May's results, d — 1). For stable matrices, the ellipses must be fully 
contained in the left half-plane (a < d), as the real part of all eigenvalues must be negative to attain stability. Accord- 
ingly, the stability criteria for fully connected matrices are: random, a\/S < d; predator-prey, (7\/S{l — 2/tt) < d; 
mixture of competition and mutualism, a\/~S{l + 2/7r) < d. To account for general connectance, C < 1, we follow 
May IT] |2l and include it under the square root, obtaining the stability criteria in Table 1 . 

We have confirmed these results by plotting the density of the eigenvalues in the complex plane (Figure 1, top): 
even for matrices of moderate size (50 species or more), the approximation is very accurate. To show the sharpness 
of the transition from high to low probability of stability, we performed extensive numerical simulations (Figure 1, 
bottom). 

In all the above cases, Mij — 0, and the expected row (column) sum is also 0. What happens if we relax 
these constraints? The most extreme case is that of mutualism (Ma = —d, whereas Mji, Mij are drawn from 
|A^(0, a^) \ with probability C and zero otherwise), in which the mean coefficient is Mij = a^plpK and the expected 
row (column) sum is i? = {S — l)CMij — d. In these matrices, we find an extreme, real eigenvalue A^ = R, 
while the remaining eigenvalues (for C = 1) are arranged in a circle centered at {—d — CMij, 0) (Figure 2, top). 
Numerical simulations suggest that the radius of this circle is approximately Kti/{2 + tt). Because for mutualism 
the stability is determined exclusively by A^ = R, the criterion becomes R < Q, which is equivalent to diagonal 
dominance lfT2]| . Therefore, in this type of matrix the fact that interactions are arranged in pairs does not influence 
stability. For the competition case, the situation is reversed (Figure 2, bottom): we find an extreme negative eigenvalue 
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Table 1: Stability criteria for different types of interactions and network structures. In all cases, the criteria hold for large, S X S 
Gaussian matrices with connectance C and diagonal coefficients —d. The elements of the matrix have mean square a^. The 
competition criterion holds for C — 1. Numerical simulations report, for a given combination of parameters, the largest S 
yielding probability of stability > 0.5 (computed using 1000 matrices). In parenthesis the analytical predictions. 




Figure 1: Top: Distribution of the eigenvalues for random, predator-prey and mixture of competition and mutualism matrices. For 
S — 250, C — 0.25 and a = 1, we plot the eigenvalues of 10 matrices (colors) with on the diagonal and off-diagonal elements 
following the random, predator prey or mixture prescriptions. The black ellipses are derived analytically in the text. Bottom: 
Corresponding stability profiles. For the random case, starting from S = 250, C — 0.5, u = 0.1 and d = 1, we systematically 
varied C (x) or a (-\-) in order obtain K = aVSC spanning [0.5, . . . , 1.0, . . . , 1.5] of the critical value for stability (indicated 
in red, 1 in the case of random matrices). The profiles were obtained computing the probability of stability out of 1000 matrices. 
The predator-prey case is as the random but with a = 0.5 and critical value n/{'K — 2). The mixture case is as the random but 
with critical value 7r/(7r + 2). 
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Figure 2: Distribution of the eigenvalues for mutualism and competition, obtained from 10 matrices ( colors), with S = 250, C = 1, 
a = 0.05 and d — 0. In both cases we observe an extreme, real eigenvalue whose value equals the row-sum in the matrices 
(circled in gray). The remaining eigenvalues, for C = 1, are contained in a circle of radius aV SCti /{2 + tt). 



= R = —{S — l)CAIij — d, and the others are contained (for C = 1), in a circle centered in (— d + CMij, 0) 
with a radius of Kit/ {2 + tt). The maximum eigenvalue is at the very right edge of the circle, so that we can derive 
a stability criterion only for C — 1: a{\/S'K/{2 + tt) + < d. For the general case in which C < 1, the 

non-extreme eigenvalues are approximately contained in an ellipse, but we have not found an exact expression for 
general C. 

Returning to the predator-prey case, we assess the effect of imposing realistic food web structure. In community 
matrices describing food webs, we expect producers to have positive columns and negative rows, with the converse for 
top predators. These variations are likely to move some eigenvalues "vertically" in the complex plane (as the large row 
and column sums have opposite signs). To test this hypothesis, we plotted the eigenvalues for predator-prey webs in 
which interactions are arranged following the cascade 1 13 1 and niche 1 14| models (Appendix). Although these models 
cannot fully reproduce food web structure l[T4l[T5l[T6l . they are widely used to simulate food webs (e.g., ifTTHTSl ). It 
is generally believed that including realistic food web structure should increase stabiHty (e.g., JT9ll20l ). 

In both cases, we find several extreme eigenvalues with large imaginary part (Figure 3). Because the eigenvalues 
of these matrices must have the same mean and variance |10| observed in the unstructured predator-prey case, we 
observe a distortion of the ellipses, yielding larger real parts than expected. Thus, cascade and niche models should 
produce networks that are less stable than their unstructured predator-prey counterpart, with the niche model having 
a larger discrepancy. These considerations are confirmed by numerical experiments (Figure 3). We conclude that, 
surprisingly, imposing realistic food web structure hampers stability. 

In the same spirit, we measured the effect realistic structures produce in mutualistic networks. Several published 
mutualistic networks are bipartite 11211 [22l [23l [24l : there are two types of nodes (e.g., plants and pollinators), and 
interactions occur exclusively between different types. Also, bipartite mutualistic networks tend to be nested ETIl : 
the interactions of the specialists are subset of those of the generalists. Nestedness is believed to beget stabiUty 
II22I [23] 1241 . We drew the eigenvalues for these two types of structure (bipartite, bipartite and nested Appendix, 
Figure 4), and contrasted the results with those obtained for the unstructured mutualistic case. We observed several 
large real eigenvalues as in the unstructured mutualistic case, but now for each positive eigenvalue we found an equally 
large negative one (Figure 4). The bipartite case yields row sums that are substantially similar to the unstructured case. 
Accordingly, we do not expect a large discrepancy in stability. However, given that nested structures will yield (on 
average), larger maximum row/column sum (associated with the generalist plants and animals), nested structures are 
inherently less stable than unstructured ones. These results are confirmed by numerical simulations (Figure 5). 
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We have so far considered how the arrangement of the interaction coefficients affects stability, we now assess 
the role of interaction strength distributions. We have extracted coefficients from normal (or half-normal) distribu- 
tions, where the majority of interactions are close to and thus "weak". These "weak interactions" are thought to 
contribute considerably to the stability and persistence of natural systems ||25]|26l . To examine their effect, we ex- 
tended our analysis to the cases in which the absolute value of each coefficient is taken either from a uniform or a 
gamma distribution, parametrized such that Afy = and Mf^ — (as in the normal case) (Appendix). Note that 
the different distributions represent different frequencies of weak interactions. Therefore, changing the distribution 
impacts MijMji, which, compared to the random case, has to be negative to increase stability and positive to de- 
press stability. When weak interactions are preponderant, MijAIji is expected to be smaller in magnitude. Therefore, 
weak interactions should increase the stability for mutualistic and competitive systems, but decrease the stability of 
predator-prey matrices and have no effect in the random case. This argument is matched by analytical and numerical 
results (Appendix), showing that, contrary to the current belief, weak interactions can be destabilizing. 

To summarize the stability properties of the various matrices, we performed numerical simulations for all the 
types of networks and with three parametrizations (Table 1). We searched for the largest S yielding a probability of 
stability > 0.5 (measured using 1000 matrices). In all cases, S increases when moving from nested mutualism to 
predator-prey. 

We have shown that arrangement into signed pairs of interactions has large impact on stability. For example, 
the random, mixture and predator-prey matrices contain basically the same coefficients, and the large difference in 
stability is driven exclusively by their arrangement (pairs with random signs, pairs with same sign and pairs with 
opposite signs, respectively). This is consistent with the fact that the variance of the eigenvalues is driven by the mean 
product between nonzero pairs MijMji. 

Mutualism, competition, and their mixture, although yielding the same pairwise mean product, and the same mean 
and variance of the eigenvalues, have very different stability properties: increasing the fraction of competitive inter- 
actions increases stability. In these cases the stability is driven by higher moments of the eigenvalue distribution. We 
conjecture that the stability properties in these matrices are thus influenced by products of three or more coefficients 
at a time. For example, the product of three competitive interactions would yield a different sign from that of three 
mutualistic interactions, potentially accounting for the difference in stability. 

Surprisingly, imposing realistic structure to the interaction networks appears to be detrimental for stability both 
in the predator-prey and the mutualistic cases. This does necessarily mean that more realistic networks should be 
less stable, as in our comparisons all the coefficients have (in absolute value) the same expectation, while this is not 
the case in natural systems (e.g., generalist species will have typically lower values for each interaction compared to 
specialists). However, we can safely assert that realistic structure alone is not contributing to stability: an increase 
in stability can be observed only if there is an interplay between the network structure and the interaction strengths. 
The fact that interaction strengths are a major determinant of stability is confirmed by our analysis showing that weak 
interactions can be either stabilizing or destabilizing depending on the type of interaction between species. 

We have shown that simply plotting the density of the eigenvalues provides qualitative insight into the stability of 
the systems. Using this graphical method, the effect of any type of network structure can be readily analyzed. 

Finally, we have found a consistent "stability hierarchy" spanning mutualism to predator-prey. Predator-prey 
interactions enable the stable coexistence of networks as large and complex as those observed empirically. 

Our results are not Umited to the stability of ecological (or biological) systems. The criteria, in fact, hold for any 
system of differential equations resting at an equilibrium point. 
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Appendix 

Construction of the community matrices 

In the main text we analyze different types of matrices. Here we detail how the matrices were constructed. In all 
cases, the parameters are: S, number of species; C, desired level of connectance; <t, standard deviation of the normal 
distribution from which coefficients are drawn; — d value of the diagonal coefficients. For each matrix, we also report 
the expected mean and variance for the eigenvalues. 
Random Matrices 

In the random case, we construct the matrices in the following way: i) for each coefficient Alij, i ^ j, we draw a 
random value X from an uniform distribution f7[0, 1]. If the value is X < C, we draw the coefficient A/^ from a 
normal A^(0, cr^). Otherwise (X > C), Mij — 0. ii) All Ala — —d. For these matrices, the eigenvalues have A = —d 
and Vai-(A) = 0. 
Predator-Prey Matrices 

i) For each coefficient Mij, j > i, we draw a random value X from U[0, 1] . ii) If the value is X < C, we draw a second 
random value Y from U[0, 1]. If this new random value isY< 0.5, we draw Mij from an half-normal distribution 
|A^(0, f7^)| and Mji from a negative half-normal — |A^(0, (t^)|, while if F > 0.5 we do the opposite, iii) If X > C, 
we assign to both coefficients, iv) All Ma — —d. For predator-prey, A = —d and Var(A) — ~{S — l)Ca'^2/TT. 
Mixture of Competition and Mutualism Matrices 

i) For each coefficient Mij, j > i, we draw a random value X from U[0, 1]. ii) If the value is X < C, we draw a 
second random value Y from U[0, 1] . If this new random value isY < 0.5, we draw Alij and AIji from an half-normal 
|A^(0, 0-^)1, while if y > 0.5 we draw both coefficients from a negative half-normal distribution, iii) If X > C, we 
assign to both coefficients, iv) All Ala = —d. For this mixture, A = — d and Var(A) = (S" — l)(7cr^2/7r. 
Mutualism Matrices 

i) For each coefficient Alij, j > i, we draw a random value X from J7[0, 1]. If the value is X < C, we draw both 
Mij and AIji from an half-normal |iV(0,CT'^)|. Otherwise, we assign to both coefficients, ii) All Ala — —d. For 
mutualism, A = -d and Var(A) = (5 - l)Ca^2/n. 
Competition Matrices 

i) For each coefficient Alij, j > z, we draw a random value X from [/[0, 1]. If the value is X < C, we draw both Alij 
and AIji from a negative half-normal — |Af(0, cr^)|. Otherwise, we assign to both coefficients, ii) All Ala = —d. 
For competition, \~—d and Var(A) — (S — l)(7(T^2/7r. 
Cascade Predator-Prey Matrices 

In the cascade model [ ISJ, species are ordered and each species has a fixed probability of preying upon the preceding 
species. The produced networks do not contain cycles, although cycles are observed in empirical networks lfT4l . In 
the cascade model, the species with highest ranking functions as a top predator, while that with the lowest ranking as 
a producer. Accordingly, the highest ranked has positive column and negative row, while the opposite is true for the 
lowest ranked. 
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The matrix construction algorithm is: i) For each coefficient Mij, j > i, we draw a random value X from U[Q, 1]. 
If the value is X < C, we draw the coefficient Mij from an half-normal |iV(0, a'^)\ and the coefficient AIji from a 
negative half-normal — |iV(0, (t'^)\- Otherwise, we assign to both coefficients, ii) All Ma = —d. For the cascade 
model, X = -d and Var(A) = -{S - l)Ca^2/TT. 
Niche Predator-Prey Matrices 

The niche model |14| allows for trophic cycles and cannibalism. The species are ordered (each one being assigned 
a "niche value", rji). A "niche radius", r^, proportional to rji, is drawn for each species along with a "niche center" 
Ci. Each species i preys upon all the species whose rjj are included in the range [ci — ri/2,Ci + ri/2]. The produced 
networks are interval (i.e., each predator preys upon consecutive species). Empirical networks, however, are not 
perfectly interval lll4l [TSl [T6l . 

To generate the matrices, we first produced an adjacency matrix A, using the niche model {Aij = 1 if i is a prey 
of j). Then we obtained a "sign matrix" S = ~A + A*. Finally, Mij is obtained multiplying Xij taken from an 
half-normal distribution |iV(0, (T^)| by Sij. The diagonal elements are set to —d. For the niche model, A = —d and 
Var(A) = -{S -l)Ca'^2/TT. 
Mutualistic Bipartite Matrices 

For the bipartite case, we divided the species in two group of equal size {S/2, when S is even). For each Mij where i 
belongs to the first group and j to the second, we draw Mij and AIji from an half-normal distribution | A^(0, cr^) | with 
probability C — 2C{S — (so that the expected connectance is matched). The diagonal is —d. For the bipartite 
mutualistic model, A = —d and Var(A) = (5 — l)C(T^2/7r. 
Mutualistic Nested Matrices 

Nestedness is a property of the incidence matrix B (typically rectangular) in which the row are the species belonging 
to the first group (e.g., plants) and the rows those in the second group (e.g., pollinators). Say that to match the desired 
connectance C in the matrix M, we want to arrange L links in B. We arrange the links in the following way: First, 
to guarantee connectedness, i) we fill the first row; ii) we fill the first column; iii) We arrange the subsequent links so 
that the matrix is perfectly nested. For example (using a squared incidence matrix), say that i? is 6 x 6: 



B = 



1,1 


1,2 


1,3 


1,4 


1,5 


1,6 


2,1 


2,2 


2,3 


2,4 


2,5 


2,6 


3,1 


3,2 


3,3 


3,4 


3,5 


3,6 


4,1 


4,2 


4,3 


4,4 


4,5 


4,6 


5,1 


5,2 


5,3 


5,4 


5,5 


5,6 


6,1 


6,2 


6,3 


6,4 


6,5 


6,6 



and suppose we want to include 17 links. First, we fill the first row (1, 1) to (1, 6), so that we placed 6 links. The next 
five links are used to fill the first column (2, 1) to (6, 1). Finally, the last six links are placed in (2, 2), (2, 3), (3, 2), 
(2, 4), (3, 3), etc. Note that the sum of the x and y coordinates for the links is growing. In fact, ordering the potential 
link by their coordinate sums, and giving precedence to those with smaller row number in case of ties, guarantees the 
maintenance of perfect nestedness. This is the filling algorithm we used in the simulations. Once we obtain B, we 
use it (along with its transpose), to determine the interactions in the matrix M. All the nonzero values of Mij are 
taken from the half-normal |A^(0, cr^)|. The diagonal is —d. For the bipartite nested mutualistic model, A = —d and 
Vai-(A) = {S ~l)Ca^2/Ti. 
Different Distributions and Weak Interactions 

The goal of this section is to investigate whether our findings are robust to changes in the distribution of interaction 
strengths and to assess the role weak interactions play for stability. In the main text, we deal with normal (or half- 
normal) distributions. Here, we consider uniform distributions and gamma distributions with different shapes. 

To measure how preponderant "weak interactions" are in a distribution of interaction strengths, we take the ex- 
pectation for the absolute value. For example, for a normal distribution, define X as a random variable taken from the 
distribution |A^(0, cr^)|. The expectation for X is E[X] — a^j2jT: w 0.798(t (as we saw in the main text). 

Before we can examine the other distributions, we have to parametrize them in such a way that Mij — and 
Mf^ = (7^. In this way, we are satisfying the conditions we stated in the main text, and we are considering distributions 
of interaction strengths with the same mean and variance. 

We start from a uniform distribution ?7[0, S\. We will sample positive coefficients from the distribution and nega- 
tive coefficients reversing the sign. Thus, the interaction strengths (for the nonzero terms) are distributed uniformly in 
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U[-0, 9]. Clearly, the distribution satisfies Mij = 0. When 6 = £7\/3, the variance M?. = ^O'^a'^ = u^. With such 
a parametrization, the expected interaction is E\X\ = (T\/3/2 ss 0.886(t: typically, interactions will be larger than 
in the normal case. What are the consequences for stability? We begin by computing r. For a predator-prey matrix 
whose elements are taken from the uniform distribution above, the expected product of the interaction strengths of 
two interacting species is MijMji = and thus r = — |. This means that the stability criterion for uniform 

predator-prey becomes a\fSC < 4d: the uniform predator-prey matrices are much more hkely to be stable than their 
normal counterpart. The uniform random and mixture cases trivially follow. For the uniform random, ay/SC < d: 
exactly as in the normal case. For the uniform mixture, a^TSC < 4(i/7, less than in the normal case. To sunmiarize, 
for the uniform distribution (where we expect stronger interactions than in the normal case), we observe an increase 
in stability for the predator-prey case and a decrease for the mixture case. The random case in unaltered. 

To further prove that weak interactions stabilize the mixture and destabilize the predator-prey cases, we analyzed 
the effect of taking the magnitude of the coefficients from a gamma distribution. The gamma distribution takes two 
parameters, k (shape), and 9 (scale). As above, we want to ensure that Mij = and Mfj = cr^. This is accomplished, 

for arbitrary k, by choosing 6 — a^ ^^^^_^^ y. For such a gamma distribution, E\X\ = For example, for 

k = 0.7, E[X] w 0.642cr (less than in the normal case), while for fc = 3, E[X] « 0.866(7 (more than in the normal 
case, and exactly as in the uniform case). Computing t for the predator-prey case, we find r = —k/{k + 1) and 
thus the stability criterion becomes aV SC < {k + l)d. This means a higher likehhood of stabihty, compared to the 
normal case, whenever k + 1 > n/{'K — 2) (approximately, k > 1.75). Thus, when we extract the magnitude of the 
coefficients from a gamma distribution, increasing k will increase the expectation, and this in tum will result in more 
stability for the predator-prey case and less for the mixture (criterion: aVSC < d{2k + l)/(fc + 1)). For the random 
case the situation is unaltered. 

Our analytical predictions are confirmed by plotting the density of the eigenvalues (Figure 6) and drawing stabihty 
profiles (Figure 7) for all the distributions described above. 
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Figure 3: Top: Distribution of the eigenvalues for cascade and niche models, with the same values used in Figure 1. In both cases 
we observe extreme, largely imaginary eigenvalues. Accordingly, the theoretical ellipse for the predator-prey case (black) does 
not contain all the remaining eigenvalues, decreasing the probability of stability. This effect is more pronounced in the niche 
case. Bottom: Stability profiles obtained using the same values as Figure 1 (unstructured predator-prey case). We report the 
unstructured predator-prey profile (blue) for comparison. Note that, contrary to the unstructured model, in the cascade and niche 
matrices increasing C or o yields different effects (separation between the two lines). This is because the row (column) sum scales 
linearly with C and a, while K scales sublinearly with C. 
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Figure 4: Distribution of the eigenvalues for the unstructured mutualism, bipartite mutualism and nested mutualism. In all cases, 
S = 250, a — 0.1, C = 0.2 and d = 0. Note that the bipartite case does produce extreme negative real eigenvalues coupled 
with positive ones, but the row sum (and thus the dominant eigenvalue) is similar to that of the unstructured mutualistic case. The 
nested matrices, in which generalist species yield { on average), larger row and column sums, display larger dominant eigenvalues. 
Thus, nestedness should produce matrices that are less likely to be stable compared to the other two cases. 
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Figure 5: Stability profile for the mutualistic cases, for S — 250, cr — 0.025, d — 1. We vary C so that the critical value 
(S — l)Ca^2/iT — R + d (x-axis) .spans [0.2, 1.2]. The critical value is reached for (S — l)C(T-^2/7r ~ 1 ( red line). Note that, 
as expected from Figure 4, nested matrices are much less stable than the other two types of matrices. 
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Figure 6: Distribution of the eigenvalues for random, predator-prey and mixture of competition and mutualism matrices (rows) 
for different distributions (top row). For S — 250, C = 0.25 and a — 1, we plot the eigenvalues of 10 matrices (colors) with 
on the diagonal. For the off-diagonal elements, the magnitudes are taken from an half-normal (first column ), gamma ( second 
and third columns) or uniform (fourth column) distributions, while the signs are assigned according to the types of matrices. The 
distributions are parametrized in a way such that Mij = and Mf- = a^. The ellipses are derived in the Appendix. 
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Figure 7: Stability profiles for the combinations of matrix type and distribution illustrated in Figure 6. For the random and mixture cases, 
stalling from S = 250, C = 0.5, a = 0.1 and d = 1, we systematically varied C (x) or a (+) in order obtain K = a\J SC spanning 
[0.5, . . . , 1.0, . . . , 1.5] of the critical value for stability (indicated in red, 1 in the case of random matrices). The profiles were obtained computing 
the probability of stability out of 1000 matrices. The predator-prey case is as the random but with cr = 0.3 for the half-normal and r(0.7, 0.917(t), 
while (T = 0.5 for the uniform and r(3, 0.289(t) cases. The adjustment of cr is necessary as it would otherwise lead to C > 1. 
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